from math import *
from numpy import *
import matplotlib.pyplot as plt
from scipy.optimize import curve_fit

print("Práctica de Dispersión de Rutherford (RU)")

N_fondo = 2   # en 600s

t_tita = [100, 100, 210, 200, 200, 200, 400, 400, 900, 900, 1200, 1200] # segundos
N_tita = [3329, 3265, 4597, 3901, 1034, 829, 330, 215, 227, 121, 99, 73]
tita = [-5, +5, -10, +10, -15, +15, -20, +20, -25, +25, -30, +30] # º

N_tita_s = [N_tita[i]/t_tita[i] for i in range(len(tita))]

N_tita_s_corregido = [float(N_tita_s[i]*2*pi*sin(abs(radians(tita[i])))) for i in range(len(tita))]   # correción espacial

# Gráfico no corregido
plt.scatter(array(tita), N_tita_s)

plt.xlabel("θ (grados)")
plt.ylabel("N(θ)")
plt.show()

# Gráfico y ajuste corregido
tita = [-15, +15, -20, +20, -25, +25, -30, +30] # º
N_tita = [1034, 829, 330, 215, 227, 121, 99, 73]
N_tita_s_corregido = [float(N_tita_s[i]*2*pi*sin(abs(radians(tita[i])))) for i in range(len(tita))]   # correción espacial

x = deg2rad(array(tita))
y = array(N_tita_s_corregido)

def ajuste(theta, A, B):
    return A/(sin((theta-B)/2)**4)

popt, pcov = curve_fit(ajuste, x, y, p0=[2000, 28])

A_fit, B_fit = popt

print("A =", A_fit)
print("B =", B_fit)

# Curva suave para representar el ajuste
theta_suave = linspace(min(x), max(x), 1000)
funcion_ajustada = ajuste(abs(theta_suave), A_fit, B_fit)

# Gráfica
plt.scatter(array(tita), y)
plt.plot(rad2deg(theta_suave), funcion_ajustada)

plt.xlabel("θ (grados)")
plt.ylabel("N(θ) corregido")
plt.show()

# Lámina Aluminio
d_AU = 2.5*10**-6
d_AL = 8*10**-6
Z_AU = 79

def Z_AL(N_AL, N_AU):
    return float(sqrt((N_AL*d_AU*Z_AU**2)/(N_AU*d_AL)))

N_AL = 1/2 *((112/1000)+(79/1000))
N_AU = 1/2 * ((253/100) + (262/100))

print("========================================================================")
print(f"Z_AL: {Z_AL(N_AL, N_AU)}")